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Manipulating expressions in many-body perturbation theory becomes unwieldily with increasing 
order of the perturbation theory. Here I derive a set of theorems for efficient simplification of such 
expressions. The derived rules are specifically designed for implementing with symbolic algebra 
ON ' tools. As an illustration, we count the numbers of Brueckner-Goldstone diagrams in the first several 

orders of many-body perturbation theory for matrix elements between two states of a mono-valent 
, system. 
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Many-body perturbation theory (MBPT) has proven to be a powerful tool in physics [l[ and quantum chemistry 0] • 
Although MBPT provides a systematic approach to solving the many-body problem, the number and complexity 
of analytical expressions becomes rapidly unwieldily with increasing order of perturbation theory. At the same 
i— i 1 time, exploring higher orders is desirable for improving accuracy of ah initio atomic-structure methods. Here a 
43 ! number of applications may benefit, ranging from atomic parity violation Q and atomic clocks 0, Q to a precision 
Oh- characterization of long-range inter-atomic potentials for ultra-cold collision studies [f| . 

S To overcome an overwhelming complexity of the MBPT in high orders one has to develop symbolic tools that auto- 
mate highly repetitive but error-prone derivation of many-body diagrams. The advantage of using symbolic algebra 
computing for these goals has been realized for a number of decades. For example, the pioneering "Schoonschip" 
program |7[ and other symbolic packages are employed for evaluating Feynman diagrams in quantum electrodynam- 
C/3 ■ ics and high-energy physics. We also note similar efforts in quantum chemistry Q (see also Ref. [j| and references 
therein). In atomic MBPT, developing symbolic tools was reported by the Notre Dame [lcj . Michigan [Tl| . and very 
recently by the Sydney [HJ and Kassel groups [T^ ]. 

In practical applications of MBPT one deals with products of strings of creation and annihilation operators. Typi- 
cally such products are evaluated with the Wick's theorem (see, e.g., discussion in Ref. [HI)- This is point of departure 
of symbolic calculations described in Refs. [l(| and The application of the Wick's theorem results in a series of 
Kronecker delta symbols. The next step in a derivation requires carrying out summation over the delta symbols. In 
a typical application the resulting terms are redundant and require additional efforts to further simplify and combine 
the expressions. A complexity of both the applying the Wick's theorem and the further simplification grows rapidly 
as the order of perturbation theory increases. 

Over the past decade, our group in Reno has developed an alternative set of symbolic tools for MBPT. The goal 
of our work was to study high orders of MBPT, e.g., fourth-order contributions to matrix elements for mono-valent 
atoms (l5l . Il6j . In our practical work, we found that the conventional approaches based on the straightforward 
applications of the Wick's theorem require prohibitively long computational times. To overcome this difficulty, I have 
1 derived a set of rules enabling efficient derivation of MBPT expressions for fermionic systems in high-orders of MBPT. 
These theorems are reported here. 

The Wick's theorem works at the level of elemental pairwise contractions of creation and annihilation operators. 
The basis idea of the present approach is to shortcut directly to the resulting expressions for a typical operations in 
MBPT, without the need to apply expensive pairwise operations. The theorems are formulated as a set of symbolic 
replacement rules, ideally suited for implementing with symbolic algebra systems. We provide an accompanying 
Mathematica package downloadable from the author's web-site [I?} • In this work we focus on mono-valent systems. 

The paper is organized as follows. In Section[lJ we review main results from the many-body perturbation theory and 
introduce notation. In Sections III Al and III CI we derive rules for multiplying second-quantized operators with atomic 
wave- functions. Similar theorems are derived for determining MBPT corrections to energies and matrix elements in 
Section [ill! Finally, as an illustration, in Section HVl we derive explicit expressions for matrix elements in several first 
orders of MBPT and count the number of resulting diagrams. 
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I. BACKGROUND AND NOTATION 
A. Second quantization, normal forms, and the Wick's theorems 

We start by recapitulating relevant notation and results from the second-quantization method as applied to fermionic 
systems. 

At the heart of the second quantization technique lies an expansion of the true many-body wave function over 
properly anti-symmetrized products of single-particle orbitals (the Slater determinants) . The machinery is simplified 
by introducing the creation (at) and the annihilation (a./.) operators satisfying the anti-commutation relations 

a l a l = - a Wj> 

djClk = —akdj, 
a 3 a k = 5 jk - a \ a J , 

djdj = , 

44 = 0. 

Applying strings of creation and annihilation operators to the vacuum state |0) builds the Slater determinants. A 
one-particle operator in the second quantization (such as an interaction with external field) reads 

Z = ^2 Z kl a \ a l ■ (!) 
kl 

A two-particle operator (such as a pair-wise Coulomb interaction between electrons) is represented as 

G = -'^2g l jkiala <! j aia k (2) 

ijkl 

= j'^gijkial a j a t a k ■ (3) 

ijkl 

Matrix elements z k i and gijki are conventionally defined on the basis of single-particle orbitals. Symmetry of the 
two-particle operator with respect to permuting electron labels leads to gijki = gjiik- By renaming summation 
indexes in Eq. (JSJ) and using anti-commutation rules, we may express G in terms of the anti-symmetrized combination 
g%jki = gijki — gijik, Eq. ©. Apparently, swapping indexes leads to the following properties, 

gijki = —gijik = —gjiki = gjiik- (4) 

Conventionally, in applications of the second quantization technique to many-electron systems, one distinguishes 
between three groups of single-particle states (orbitals): core, virtual (excited) and valence orbitals. The core orbitals 
are occupied and form the quasi-vacuum state |0 C ). Virtual orbitals are unoccupied in |0 C ). We will treat valence 
orbitals as a part of the set of virtual orbitals. We follow a convention of Ref. [lj] and label core orbitals as a,b . . ., 
excited (virtual) orbitals as m, n, and valence orbitals as v, w. Indexes i,j,k,l run over both core and virtual 
orbitals. For example, in the independent-particle-approximation, a state of a mono-valent atom may be represented 
as aj|0 c ), where the quasi-vacuum state |0 C ) = (rLecore !) |0) represents a closed-shell atomic core. 

Further, we review several results related to the normal form of operator products, : ••■ :. The operators are 
rearranged so that a COIC and a\ h . t appear to the left of al ore and a v irt ■ When acted on the |0 C ) quasi-vacuum state, 
most of the strings of operators in the normal form produce vanishing result. 

One of the central results related to the normal forms is the Wick 's expansion into normal products 

A=: A: + :A:, 

: A : being a sum of normal ordered terms obtained by making all possible single, double, triple ... contractions within 
A. Contractions between two creation/annihilation operators x and y are defined as xy = xy— : xy :. The overall 
sign of : A : is (— l) p , p being a number of permutations to bring A into the normal form. The same rule holds for 
terms in : A : - we count permutations necessary to bring the operators being contracted next to each other and also 
the permutations needed to arrange the resulting term in normal form. The only nonvanishing contractions are 



aradn — ^mn, CL a CL^ — O a }). 
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All contractions between the core and excited (including valence) orbitals vanish. 

Using the Wick's expansion, one can rewrite a single-particle operator Eq.fl]), as a sum of zero-body (scalar) 
and one-body contributions [14| 

Z = Z + Zi , (5) 

Zq = ^ Z aa , 



Zij : a\a,j : 



Similarly, any two-particle operator G, Eq.©, may be represented as a sum of zero-body G^- ', one-body G^ and 
two-body terms, 



G = G + Gi + G 2 , (6) 
1 
2 



G*0 — ^ ^ ffabafc 



'J ' ' 



Gi = J2 ( J2 ~ bia ^ ) : a i a 

ij \ a / 

G2 = - ^2 hjki ■ aja^aiak ■ ■ 

ijkl 

Technically, the MBPT formalism requires multiplying second-quantized operators. Simplification of the resulting 
expressions is greatly aided by the Wick 's theorem for operator products. For two strings of operators in the normal 
form : A : and : B :, the theorem states 

: A : : B :=: AB : + : AB : . (7) 

Here : AB : represents the sum of the normal-ordered terms with all possible contractions between the operators in 
A and those in B. 



B. MBPT for mono-valent atoms 

To gain insight into a general structure of MBPT expressions, here we briefly reiterate MBPT formalism [l5| for 
atoms with a single valence electron outside a closed-shell core. The lowest-order valence wavefunction is simply 
I'i'i ''} = aj|0 c ), where v is a valence orbital. The perturbation expansion is built in powers of the residual interaction 
Vi defined as a difference between the full Coulomb interaction between the electrons and the model potential used 
to generate the single-particle basis. The n th -order correction to the valence wavefunction may be expressed as 

\yW) = -R v \QV I \Vi n -V)} , (8) 

J linked 

where R v is a resolvent operator modified to include so-called "folded" diagrams [l5| . projection operator Q = 
1 — |5'i ^)(4'i ' ) |, and only linked diagrams are to be kept. For mono-valent systems a convenient starting point 
is a single-particle basis generated in the frozen-core (V N ~ 1 ) Hartree-Fock (HF) approximation In this approxi- 
mation, the residual interaction is simplified to a two-body part, G2, of the Coulomb interaction and the number of 
MBPT diagrams is substantially reduced [ToL [T^j. 

The recursion relation, Eq. ([8]). systematically solves the many-body problem, as we may generate corrections to 
the wave-function at any given order of the perturbation theory. With such calculated corrections to wave-functions 
of two valence states w and v, the n th -order contribution to matrix elements of an operator Z may be determined as 

n-l 

v(n) X^/^(n-fe-l)|^||r,(fe)\ 1 7(n) (q\ 

^wv — / , \^w \^\^v /val, conn i wv, norm ' \ u I 

k=0 

Here Z^v, norm is a normalization correction arising due to an intermediate normalization scheme employed in deriva- 
tion of Eq. ([5]). Subscript "val, conn" indicates that only connected diagrams involving excitations from valence 
orbitals are included in the expansion. 
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Generic contribution to wavefunction 



Now we would like to introduce short-hand notation for strings of creation (a k ) and annihilation (a^) operators in 
the Fermi statistics. String of x operators 

E i = 44 • • • 4 

combines creation operators for excited orbitals and symbol a ranges over the set 1, 2 ... x. Similarly, 

Cp = a\a2 ■ ■ ■ a y 

represents a string of y annihilation operators for core orbitals, with symbol (3 spanning the indexes 1,2...?/. Finally, 
V' is either aj or 1 depending on the presence of the valence creation operator. 

On very general grounds a generic piece of atomic wave-function for a mono-valent atom may be represented as 



|$)= Y, L[(l,2... x) a , (1,2... y) ] EtC Vi\Q c ), 

M,{/3} 



(10) 



where ^ 



{a} 



Ei-'-E^andL (1,2. (1,2. ..y). 



is a c-number object which depends on the indexes. 



As an illustration, a doubly-excited core state of a mono-valent atom may read ^ mnab Pmnaba\n a n a a a b<A)\^c) 



Apparently, x 



2, El 



a* at 



a a ab, V* = aj and the groups of indexes are (1,2... x) a = (m, n) and 



(1,2... y)o = (a, b). In summations, the symbolic index a would run over symbols m and n (which range over all 
excited orbitals). Similarly, the symbolic index j3 assumes symbolic values a and b, i.e., the labels of the core orbitals. 

The Eq. IjlOp will be the central object in our derivations presented below. We will act on this "generic piece of 
wave-function" with one- and two-body operators and also bring the resulting expressions into the very same form 
ofEq. m. 



II. SIMPLIFICATION THEOREMS 



We would like to follow the MBPT prescription ([8]) and compute the wave- function in an arbitrary order. To this 
end we will derive rules for simplifying products of one- and two-particle operators with the generic wave-function 
<!>}, Eq. (|10[) . The wave- function |$) is in the normal form with respect to the quasi- vacuum state. Therefore, we 
take the second-quantized operators expanded in the normal forms, Eqs. (|5l6[) . and apply the Wick's theorem ([7]). 
Apparently, the zero-body terms Zq and Go (c-numbers) do not produce non-trivial results and in the derivation 
below we focus on contractions of one- and two-body operators with strings of operators entering |$). Finally, we 
bring the resulting chain of operators in the same standardized sequence of operators E^ a Cp V* , as in (fTU]) ; of course 
the number of operators in each group may differ from the starting numbers of operators in |$). 



A. Product with a one-body operator 

Here we focus on acting with the operator 

Z x =Y z (hj) ■ 4°j : ( n ) 

on a generic wave-function |$), Eq. (|10p . According to the Wick's theorem, we may have 0, 1, and 2 contractions 
between Z\ and operators entering |$). There are only 6 distinct possibilities classified by the number and type of 
contractions 

{Zx |$)} = {Z l |$>} 

+ {^ii$>}i e + {^ii < mi c + {^ii < i>}}i, 

+ {^l| < f)}le.l C + {^l| ( f ) )}l,,lc • 

The first term corresponds to no contractions. The term le results from contracting one excited orbital from Z\ 
and one orbital from the string E\ t of the operators entering |<I>). Subscript lv, lc labels double contraction: one 
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contraction lv involves an operator from the valence string and the other contraction lc involves an operator from 
the core string Cp. The labeling scheme for other contributions follows from these examples. 

Notice, that a straight-forward application of the Wick theorem ([7]) is inefficient. Operators resulting from expanding 
Eq.([7]) ultimately act on |0 C ) and many terms in the Wick's expansion will produce zero result. Indeed, we may write 
explicitly 

Z\ = Y~^ z(m, n) a^a n + Y~^ z(a, m) a\a m — YJ^ z(a, b) a b a\ L + YJ^ z(m, a) aj„a a . 

rein ma ab ma 

As an example, consider term with no contractions. Only the last contribution from the above expansion will 
contribute, because the first two terms annihilate an unfilled orbital in |0 C ) and the third term promotes an electron 
into already occupied core orbital. Based on this discussion, we shortcut the application of the Wick theorem and 
limit ourself to a much smaller subset of terms. 



No contractions 



In this case, the relevant part of the operator Z\ contains one core annihilation operator and one creation operator 
involving excited orbital. Then we may simply move the operators from Z\ to the ends of the excited E* a and core Cp 
operator strings in the wave-function. The additional core and excited orbital indexes are absorbed in the summation. 

{Z x |$)} = YJ L° [(1, 2 . . .x,x + l) a , (1, 2 . . . y, y + 1)J (Cpa y+1 ) |0 C ) , 



where 



L° 



(l,2...x,x + l) a ,(l,2...y,y + l) = (-!)»* {x + 1, y + 1) L (1, 2 . . .x) a , (1, 2 . . .y) 



Moving a pair operators together does not produce a phase. The phase factor appears because the excited orbital 
operator was additionally moved through (anti-commuted with) y operators in the Cp string. 



Single contraction with the creation operator in the E' a group: It 



There are x such contractions ( aja' a = Sj a ); we contract the operators in turn. By contracting with a\ we obtain a 
string of operators a\- ■ ■ a\, and we bring the a\ to the end of this string thus acquiring a phase (— l) x 1 . By renaming 
the dummy summation indexes we may bring the resulting sequence of operators to the same form a\- ■ ■ a\a\. For 
example, as a result of contracting with aj, we obtain 



Now we rename (p 



Finally, 



^ (-1) M (-1) 



t t 
a a' 



1"2 



t t t„t 

2 fj,-l a H-l a x a i ■ 



1) and bring the resulting string into the form a 2 ■ • • a\,a\ 
S ji M)^ 1 (-i)^ 1 4 4 • • • a i-i a t-i a l4 = 

in (-lf-^-ir 1 (-i)" ai-'l-A'U'te 

c ( 1\X-1 t t t t t t 



{Zm} le = E L le [(2...x,x + l) a ,(l...y) ] (4-"44+i) CpVl\0 B ), 



with 



L 



lc 



(2...x,x + l) a ,(l,2...y) f3 \ = (-If 1 Y,z(x+l,l e )A le L\(l e ,2...x) a ,(l,...y) 



6 



where the operator Ai^ anti-symmetrizes L over the first excited index, i.e., 



AiL 



(l e , 2... =L (l e ,2...x) a ,Q f} -L (2,le,3,...s) a ,() fl — L (3, 2, l e , . . .x) a , Q p 



or 



As an example, 



A/(l,2..,z) = /(l,2..,z)-^/(l,2.., M ,...z) 



A~L e Pmn,ab — Pmn,ab Pnm.ab ■ 



(12) 
(13) 



Computationally, in symbolic algebra implementations, the symbol replacement operations in Eq. (|13[) are efficient. 
We note that Ai c L — x L. 



3. Single contraction with the annihilation operator in the Cp group: 1c 



The derivation is similar to the previous case; a\ap = Sip. There is a phase factor (~l) x+y due to the transfer of <Zj 
through E^ a Cp and an additional factor (— l) x due to moving the a\ to the beginning of the Cp group, resulting in 
the total phase of (— l) y . The result reads 

{Z|$)} lc = L lc [(l...x) a ,(2...y,y+l) f3 \Ei (a 2 ■ ■ ■ a y a y+1 ) 7+ |0 C ) , 



(1 . . . x) a , (2 . . . y, y + 1) J = (-If z (!c V + 1) A c £ [(1 . . . z) Q , (l c , 2, . . . y) 



The operator Ai c anti-symmetrizes L over the first core index and is defined similarly to Ai^ 



4- Single contraction with the valence creation operator: lv 



Contraction ajOv = Sj v . Phase (—l) x+v arises due to bringing dj to the end of E^ a Cp and the phase (— \) x due to 
transferring a\ to the end of E^ a . The result reads 

{Z\*)}iv= E L lv [(l...x,x + l) a ,a.--y) p } (a\---at +1 ) Cp |0 C >, 

M,{/3} 



with 



L 



(1. ..x,x + l) a , (1,2. ..y)g = (-l) y z(x + l,v) L (1, . . . x) a , (1, . . . y) 



5. Double contractions: one excited and one core operators (le,lc) 

Here the number of operators in both core and excited orbitals strings is reduced by one. The derivation is similar 
to the le case. 

{Z\mx e ,ic= E £ le ' lc [(2,...z) Q ,(2,...y)J (4 •••4) (a 2 ---a y ) V* |0 C >, 

{"M/3} 



L 



le.lc 



(2, . . . aO a , (2, . . . y^ = (-lp 1 * (lc le) (l e . . . i) a , (lc, • ■ ■ y) 
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6. Double contractions : valence and one core operators (lv,lc) 



{^>Wc= E i 1 "' k [(l,..i,(2,...!/)j4(«r-«^ |0c>, 

{"M/3} 



- Xv.Xc 



(1, . . . a:) a , (2, . . . y) J = (-1)» ]T z (l c , «) ^i e L [(1 . . . a:) Q , (l c , 2, . . . y) 



B. Example: zeroth-order Hamiltonian 



To illustrate the derived simplification rules consider a zeroth-order Hamiltonian in the second quantization, 

H = ^e l : a\a,i : . 

i 

Apparently, this is a special case of a one-body operator, Eq. (fTTj) . with z(i,j) = £j Sij. While the result of applying 
Hq to <E> is trivial, arriving at it via application of the derived rules is instructive. We would like to show that 

flol*) = (SaEa-S/SE/J + ( 14 ) 

Here S v — if there no valence operator present in |$) and S v = 1 otherwise. E a is a sum over all indexes in the 
string and ~Ep is a sum over core indexes in the Cp group. 

Because the matrix element z(i,j) = £j (5^ is diagonal, the only contraction classes which contribute are le, lc, and 
lv. For example, consider the case of le. We deal with the object 

(2 . . .x,x + l) a , (1,2 . . . y)J - (-l) 3 ^ 1 ^ e x+1 5 le , x+1 A u L \(l e , 2 . . . x) a , (1, . . .y) 



= {-l) x 1 e x+1 A x+1 L (x + l,2...ar) a ,(l,...j/)« 



Then 



{^o i$)} le = (-if 1 E E ^+1^+1^ [(^+1.2... x) a , (i, . . . (4 • • • 44 +1 ) yt | 0c ) . 

{«} {0} 

Further rename x + 1 — > 1 and place the o{ at the beginning of the string: 

{H |$)} le = E E [(le, 2 . . . x) Q , (1, . . . y) J (a f le • • • 4) C> ^ |0 C ) . 

On expanding the partial anti-symmetrization we encounter terms 

E E - £ ^ L [o*. 2. ie, ..») Q ,(!,■■■ y) J (4 e • • • 4) cp y f io c ) . 

{«} 1/3} 



Rename 1„ 



Therefore, 



/x and swap the operators 

J2 E +e » L [o-» 2 ' ••> * - x ) a > (*> • • ■ f) J ( a l • • • a *) c /5 Ft 

{«} {/3} 



{H \$)} u = (S a e„)L[(l,2...x) a ,(l,2...i/) /S ] E+C^yt| 0c ). 

Here E Q is a sum over symbols. Similarly, 

{JTo |$>} lc = - £ (Z p e p )L[(l 7 2...x) a ,(l,2...y) f3 \ Cp |0 C ) , 

andjiJo — Et>| < £ > )- All the remaining contractions vanish because they involve matrix elements between core 

and excited states. Finally, by adding the derived terms we arrive at the well-known formula |]2 
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C. Contractions with two-body operator 

Now we consider products of the two-body part of a two-particle operator G with our generic piece of the many-body 
wave function, Eq.©. We will use 



G 2 = j^2g(i,j,k,l) : at a] 



a t a k : . 



ijkl 



The derivation is similar to the one-body case of the preceding Section. Here, however, the maximum number of 
possible contractions is four and there are 15 distinct cases, enumerated below. 



1. No contractions 



{G 2 |$)} = L ° 

M,{/3} 



(1, . . . x + 2) Q , (1, . . . y + 2) p ] (Eftat +1 at +2 ) q (Ca y+1 a y+2 ) V* |0 C > 



L° [(!,... x + 2), (1,... y + 2)] = -g(x + l,x + 2,y + 2,y + 1) L [(!,... x),(l,...y)] 



2. Single le 



{G 2 |<f)} le - ^ 
{*}M 



(2,...x + 2) Q ,(l,...y+l) /3 ] (4---4+i4 +2 ) a (C'a y+1 ) V* \0 C ), 



L le [(2,...x + 2),(l,...y + l)] = (_i)^- 1 \Y.~9^ x + 1 > x + Me,2/ + L [(1 B) . . . x) , (1, . . . y)\ 



The partial anti-symmetrization operator A\ c is given by Eq. (|13D . 



3. Single lc 



{G 2 |$)} lc = £ i lc [(l,...^ + 2) Q ,(2,...y + 2) /3 ] (i<&4 +1 ) (<z 2 . . . a y+1 a y+2 ) |0 C 



(l,...z + 2) Q ,(2,...y + 2)J =-^3(x + l,l c ,y + 2,y + l)^l 1(! L[(l ) ...ar),(l c ,...y)] 



^. Single lv 



(l,...a; + 2) a ,(l,...y+l)J (^4+i4+ 2 ) q (G^+i)^ |0 C ), 



I 1 



(l,...x + 2) Q ,(l,...y+l), 



= -g(x + l,x + 2,v,y+i)L[(l,...x) ) (l t ...y)] 
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5. Double 2e 



{G2 |$>} 2e = £ L 2e [...] (4 •• -4+14+2) C-ZJ^lOc), 



L 2e [...] = -^j( 3 ; + l,x + 2, l e> 2 e ) ^l lc2c L [(l ei 2 e . . . x) , (1, . . .y)] . 

1 C ,2 C 

The partial anti-symmetrization operator is defined as 
A h2 f (l,2..x) = /-^/(l I 2,..7 ?I ..,ar) 77Ml - £)/ (1,2,.. 77, .., z)„ M 2 + ]T / (1, 2, ..!/, ..»,, .., ar)^^^ (15) 

r)>2 r;>2 ))>^>2 

(this anti-symmetrization produces x(a; — l)/2 terms). 

As an example, Ai c2( ,p m n = Pmn and -4i E 2 E Pmnr = Pmnr — Prnm- Alternatively, this definition may be rewritten as 

A 1>2 f (1, 2..x) =Aif + A 2 f-f + 2 /(2, 1, ...x) + Y, f 2 < - a0*~i,„~2 ■ 

?7>j/>3 

6. DowfcZe ,3c 



{G 2 |$)} 2c = E L 2c [...}(E^) a (a 3 ...a y+1 a y+2 ) fl V^\0 c ), 
L2C [■■■] = ^ E 5 (2- 1c, y + 2, y + 1) A lc2 L [(1, . . . x) , (l c , 2 C . . .y)} , 

lc,2 c 

where the partial anti-symmetrization operator Ai c2c is defined similarly to A\ c2c . 

7. Double le,lc 



{G 2 |$)} le , lc = £ L le - lc [...] (4---4+i) (a 2 ...a, +1 )^t|0 c ), 
L le < lc [...] = (-lf+y +1 Y~9(x + 1, lc le, y + 1) Xi e A le L [(l e , . . .a,) , (l c , . . . y)} , 

lc,lc 

Here we deal with two successive applications of the partial anti-symmetrization operators introduced in Section Hi Al 
Eq.lHSJ. 

8. Double lc,lv 



{G 2 |*>} 1Vilc = £ L lv ' lc l--} (^4+i) (a 2 ...a v+ i) /J |O c ), 
L lv - lc [...]= J2 9 (x + 1, l c , «, y + 1) -4i c £ [(1, • • ■ x) , (l c , . . . y)] . 



9. Double le,lv 

{G2 m Utlv = E L le ' lv l---} (4---4+i4 +2 ) c p \o c ), 

L le ' lv [...] = (-If +V Y,9 (x + l,x + 2, v, l e ) ^ le L [(l e , ...*), (1, ... y)] 

{G 2 |$)} 2e , lc = E i 2e ' lc [---] (4---4+i) (a 2 ...a v )^t|o c ) ) 

L 2e < lc [...]= £ S(x + l,lc,le,2e)^i c ^i a ,2.£[(l,...a;),(lc,...l/)]. 

iJ. JHpZe le,2c 

{G 2 |$)} le , 2c = E i le>2c [---] (4---4) (aa-.-^+O^tlOc), 

{«},{/?} 

i l6 ' 2C [.-.] = (-1) X+!,+1 E 3 ( 2 « le. 2/ + 1) ^le^lc,2 c i [(1, • • • *) , (lc, • • 

i2. IHpZe iv ; ,Sc 



{G 2 |$}} 1V;2c = E ^"^ [• • -1 ( Et ) a ^ ■ ■ ■ \°oh 

{"M/3} 

L lv,2c [...]= E 5 (2c lc, «, y + 1) ^i e , 2c £ [(1, • • • x) , (l c , . . . y)} . 

1 C ,2 C 

i3. IHpZe le,lv,lc 
{G 2 |<f)} 1Mc , le = E i lo,lc,le [---] (4---4+i) (a 2 ...%) /3 |O c ), 

i^' ic ' ie [...] = (~ir +v+1 E-9(^+ 1. ic ic, «) ^i e ^i.£ [(1, . . . x) , (i c , . . . 
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14- Quadruple 2e,2c 



{G 2 |$>} 2e , 2c = E i 2e - 2c [. 

{<*},{/)} 



(4 •■•4) ( a 3--- a v) v f \°c), 



L 2e > 2c [...]= £ g (2 C , l c , l e , 2 e ) A e ,2^i c , 2c i [(1, ■ • • z) , (l c , . . . y)} . 



15. Quadruple le,lv,2c 



{G 2 |*)} lt)lle)2c = E (4... 4) (aa...^ |0c), 

L 2e ' 2c [...] = E 9 ( 2 - lc le, «) ^i e A c ,2 c £ [(1, • • ■ x) , (l c , ...y)] 

le,lc,2c 

D. Additional remarks 



Note that the introduced partial anti-symmetrization operators, Eq. ()13|15|) are subsumed into the following general 
definition of a partial anti-symmetrization operator A [{1, £} 7 , {1, fj,}s] over a subset {1,£} of excited indexes and a 
subset {1,^} of core indexes 

A [{1, 07, {1, M>«] i [(1, 2 • • ■ x) a , (1, 2 . . . y)^ = 
^^{r[(U...i]F[(l,2...JL[(l,2...i,(l,2... !/ )J} 

Here the summation is over all possible renaming of indexes in the groups. Individual phases of the terms are 
determined by functions F which we use to denote the conventionally-defined signature of the resulting permutation 
of indexes. 

The derived rules may be presented in a more symmetric form by noticing that no matter how simple or complicated 
the dependence on the indexes inside the object L is, by systematically swapping the dummy summation indexes the 
generic piece of MBPT wave- function, Eq. (JXOJ) , may be rewritten as 

l $ ) = i E Z[(l,2... a O a ,(l,2...y)J ^C^t| 0c)) 



x\y< 



where the object L was obtained by a complete anti-symmetrization of L inside of the groups of excited and core 
indexes. In other words, in the derived theorems one could simply replace 

£[...]- J-Z[...] 
x\y\ 

and unfold the partial anti-symmetrization operators. While the resulting expressions may be more aesthetically 
appealing, we did not find any particular advantage of using them in practical calculations. 



III. OBSERVABLES 



Application of the derived rules allows us to find many-body correction to atomic wave-function in an arbitrary 
order of MBPT via Eq. ([5]). Below we focus on an efficient symbolic evaluation of expressions for MBPT correction 
to energies and matrix elements. 
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A. Corrections to energy 

The correlation correction to energy in the n-th order of MBPT may be found with the n — 1-th order correction 
to the wave-function 

5E^ = (yW\Vi\V( n -V), (16) 
or in the frozen-core approximation for mono-valent atoms, 

5EM = <0ck {g 2 |^"- x ))} . (17) 

Now we focus on the object in the curly brackets, \4>) = G^l^i™ ). We may use the results of Section lll Cl and derive 
a multitude of the terms on the r.h.s. Ultimately determination of the energy correction is simplified by noticing 
that only a small number of terms would remain after forming the required product (O c |a„|0). The non- vanishing 
contributions arise from generic pieces 

{G 2 |*("- 1 ))}^=L c [()()]at| 0c ) ! 

where L c does not depend on the valence orbital, and 

{G 2 \¥ n ~V)} =££4(l e )0]alj0 c >, 
e i e 

where L v necessarily depends on the valence index. Accordingly, the corrections to the energy may be separated into 
the core and valence parts, 

with 5Ec^ — L c and SE^ = L v [(«)()]. This solves the problem of finding the energy correction. 



B. Matrix elements of a one-body operator Z 



Suppose we derived the MBPT corrections to wave- functions of two valence states w and v. We would like to 
compute the matrix element of some one-particle operator Z. To this end we need to use the formula |9j) and compute 
terms 



We start by introducing an intermediate state 



l^^z^). 



Then our task will be accomplished by forming the product k ^\if) v )' 

\4> v ) may be computed using rules of Section [II Al It will contain a linear combination of various "generic contri- 
butions"; we focus on a scalar product of generic contributions to \ip v ) 

|$,)= ]T L [(1,2. (1,2...^] EfcCpVi \0 C ) 



and to 



(n-fc-l) 



(*™ | = J2 K[(l,2... x') a , ,(1,2... y') ,\ (0 C | WC\, E a , . 
{«'},{/9'} 

We adopt a convention that in the groups of operators C^, and E a i the enumeration of operator indexes goes from 
right to left, i. e., Cjj, = a) , ■ ■ -a\,a\, and E a > = a x > ■ ■ -a^ay. 
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To obtain a nonzero value for the scalar product ($> w \$> v ) we have to perform a maximum number of contractions. 
Moreover the numbers of operators must be related as x + $ v = x' + 5 W and y = y'. 

We contract the core orbitals first. Bringing the two groups of operators together introduces a phase factor of 
(-l)( x - x ') y and we obtain 

(* w |*„> = (-l)( x - x > E J2(0c\WE a ,Elvi\0 c ) 

{a'}, {a} {/3} 

L [(1, 2 ...*)„, (1, 2 .. . y) J „4 {/3} X [(1,2... x\, ,(1,2... y) 

Here -4{/3} denotes a complete anti-symmetrization over all indexes in the symbolic set {/?}. The result may be 
proven by noticing that the complete contraction between the products Cp, Cp may be expressed as 



CtCp = det 



$w 5w ■ ■ ■ Si x i 

$21' $22' ■ ■ ■ $2x' 
$xl' $x2' ■ ■ ■ $xx' 



The result of computing the remaining product (0 C | WEa'E^ |0 C ) depends on whether we have the valence 
operators present inside the W and groups; we treat three separate cases below. 

1. Both valence operators are absent 

In this case, both objects, L and K, necessarily depend on their respective valence indexes and we emphasize these 
dependencies by super-scripts L v and K w . Also x = x' and (similar to full contraction for core orbitals) 



(*™|^>vai = E E U [t 1 ' 2 ...*)„, (1, 2 ... 

{«} 1/3} 

Such a combination produces x\y\ diagrams. 



(l,2...x) a ,(l,2...y) fj 



2. Only one valence operator is present 



Consider first the case when only the operator a\ is present, then the object K depends on the index w; we denote 
this dependence as K w . Also x' = x + 1. By appending a\ at the end of the E^ string, we reduce the treatment of 
contractions to the preceding case. 

(*™|*„> - (-1) V E E^l E <*' E « L [(L 2 ■ • • x )a > (1 • • • f) J [(!'. 2 ' • • • ,(!■■■ ^ 

{a'},{a}{/3} 



(*™|*„>„ = (-If EE L 2 ••••<•:•„ J A"} 



If only a w is present, then 

(* w \* v ) w = (-1)" E E K [(!'' 2 ' • • • *')« , (1 • • • 2/); 

{«'} {/3} 



-4{a}-4{/3}£" 



(1,2. ..£,<,(!...</), 



(l',2'...x',w) a ,(l...y). 



3. Both valence operators are present 



If both valence operators are present, i.e., W — a w and W = aj, then x = x' . The product (0 C | a w E a /E^ a\ |0 C ) 
may be broken into two contributions. 
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(i) Contraction between a w and leads to a core contribution 
i^ w \^)T x ' = <U E E L t 1 ' •••*)„> (L •• • Vh 

M {/?} 



(l,...ir) Q ,(l,...2,), 



(18) 



Notice that this contribution vanishes for w ^ v and, moreover, it does not vanish only in a very special case of true 
scalar operator Z. 

(ii) Simultaneous contractions between a w and an operator in E^, aj, with an operator in E a / and residual contrac- 
tions lead to (ylxl) x contributions 



£ — >w 



(l,2...z) Q ,(l,...y) 



Here S represents a summation over all possible simultaneous replacement of index £ by w in L and by v in if. 



IV. COUNTING DIAGRAMS AND SUMMARY 



The present paper provides symbolic prescriptions to aid in deriving many-body diagrams for mono-valent systems 
in an arbitrary order of MBPT. Based on the derived rules, the author has developed a Mathematica package, which 
is made available through the author's website [TiTl . The rules have been tested by recovering known results for matrix 
elements and energies through the third order [lCj . Notice that the package and the theorems have been already used 
in deriving matrix elements through the fourth order of MBPT [151 1 1 61 ] . 

The derived rules are certainly not as mnemonically elegant as the original Wick's theorem; this is a reflection of 
the fact that the Wick theorem is formulated in terms of the pair-wise contractions between the operators, while our 
rules provide an explicit (yet general) answer that avoids a multitude of elemental pair-wise contractions. In this 
sense, the derived theorems may be called "post- Wick" theorems. 

As an illustration, we derive the MBPT diagrams for matrix elements of a one-body operator between two distinct 
states w and v of a mono-valent atom. We run the package and obtain analytical expressions for the diagrams. The 
expressions are too lengthy to be presented here (see, however, a partial list of fourth-order diagrams in Ref. [HI]); 
instead we count the number of the resulting diagrams. The results are compiled in Table |U The counts do not 
include the normalization term of Eq. ([9]) and the core contribution (fT8| . 



TABLE I: Complexity of MBPT for mono-valent systems in the frozen-core approximation. We list numbers of diagrams in 
the n th order of MBPT for wave-functions and matrix elements of a one-particle operator. There are two counts for Zwv in the 
format nx/n^. Ui is a number of terms in a maximally simplified expression (where the Coulomb integrals gijki were combined 
into the anti-symmetric combinations gijki)- n% is the full number of the conventional Brueckner-Goldstone diagrams including 
exchange diagrams. 



Order n \fj n) ) Zj n J 

1 

1 2 1/1 

2 20 2/4 

3 561 30/84 

4 26700 552/3072 



As discussed, due to our short-cutting the expensive elemental pair- wise contractions, the advantage of the present 
formulation becomes most pronounced in high orders of perturbation theory, as both the number of resulting terms 
rapidly grows (see Table HJ and the length of operator strings entering |\I>( n )} increases with n. In addition, from a 
practical standpoint, in symbolic algebra implementations, the next step involves combining similar terms by pattern 
matching. This is a computationally expensive search operation. By contrast to the directly applied Wick's expansion 
method, our derived theorems already yield simplified results further speeding up the symbolic evaluations. 
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